The objective of this notebook is to perform a basic quality control of the mapping performed with Cellranger. Cell Ranger ARC analyzes data generated by the Chromium Single Cell Multiome ATAC + Gene Expression. Should note here, we will talk about features that can refer to a gene (declared from the hg38 reference transcriptome) or a peak (determined by the Peak Calling algorithm provided by Cell Ranger).
The expected values (EV) included in the notebook are obtained from this technical_note
library(ggplot2)
library(ggpubr)
library(plyr)
library(reshape2)
library(data.table)
library(knitr)
library(kableExtra)
library(ggrepel)
library(plotly)
path_to_multiome <- here::here("multiome/results/tables/cellranger_mapping/cellranger_mapping_metrics_multiome.csv")
path_to_proj_metadata <- here::here("multiome/1-cellranger_mapping/data/tonsil_atlas_metadata_multiome.csv")
barplot_data <- function(multiome_melted){
multiome_melted$value <- round(multiome_melted$value, 2)
ggbarplot(multiome_melted,
x = "library_name",
y = "value",
fill = "gray70",
ggtheme = theme_pubr(x.text.angle = 45,legend = "none")) +
facet_wrap(ncol = 1,scales = "free_y", ~variable)
}
correlation_plot <- function(multiome_df, variable1, variable2)
ggscatter(qc_samples, x = variable1, y = variable2,
add = "reg.line", # Add regressin line
add.params = list(color = "blue", fill = "lightgray"), # Customize reg. line
conf.int = TRUE # Add confidence interval
) + stat_cor(method = "pearson")
qc_samples <- read.csv(path_to_multiome)
metadata <- read.csv(path_to_proj_metadata)
qc_samples$library_name <- revalue(
qc_samples$Sample.ID,
c(
"ulx1v6sz_8a2nvf1c" = "BCLL_8_T_1",
"wdp0p728_jf6w68km" = "BCLL_8_T_2",
"co7dzuup_xuczw9vc" = "BCLL_9_T_1",
"qmzb59ew_t11l8dzm" = "BCLL_9_T_2",
"pd9avu0k_kf9ft6kk" = "BCLL_14_T_1",
"vuuqir4h_wfkyb5v8" = "BCLL_14_T_2",
"admae8w2_89i88tvv" = "BCLL_15_T_1",
"sr20954q_yiuuoxng" = "BCLL_15_T_2",
"kmbfo1ab_ie02b4ny" = "BCLL_2_T_1",
"ryh4el3i_biv0w7ca" = "BCLL_2_T_2",
"bs2e7lr7_mdfwypvz" = "BCLL_2_T_3"
)
)
knitr::kable(qc_samples, row.names = TRUE) %>%
kable_styling(
bootstrap_options = c("striped"),
full_width = T,
font_size = 12
) %>%
scroll_box(height = "300px")
| Sample.ID | Pipeline.version | Genome | Estimated.number.of.cells | ATAC.Confidently.mapped.read.pairs | ATAC.Fraction.of.genome.in.peaks | ATAC.Fraction.of.high.quality.fragments.in.cells | ATAC.Fraction.of.high.quality.fragments.overlapping.TSS | ATAC.Fraction.of.high.quality.fragments.overlapping.peaks | ATAC.Fraction.of.transposition.events.in.peaks.in.cells | ATAC.Mean.raw.read.pairs.per.cell | ATAC.Median.high.quality.fragments.per.cell | ATAC.Non.nuclear.read.pairs | ATAC.Number.of.peaks | ATAC.Percent.duplicates | ATAC.Q30.bases.in.barcode | ATAC.Q30.bases.in.read.1 | ATAC.Q30.bases.in.read.2 | ATAC.Q30.bases.in.sample.index.i1 | ATAC.Sequenced.read.pairs | ATAC.TSS.enrichment.score | ATAC.Unmapped.read.pairs | ATAC.Valid.barcodes | Feature.linkages.detected | GEX.Fraction.of.transcriptomic.reads.in.cells | GEX.Mean.raw.reads.per.cell | GEX.Median.UMI.counts.per.cell | GEX.Median.genes.per.cell | GEX.Percent.duplicates | GEX.Q30.bases.in.UMI | GEX.Q30.bases.in.barcode | GEX.Q30.bases.in.read.2 | GEX.Reads.mapped.antisense.to.gene | GEX.Reads.mapped.confidently.to.exonic.regions | GEX.Reads.mapped.confidently.to.genome | GEX.Reads.mapped.confidently.to.intergenic.regions | GEX.Reads.mapped.confidently.to.intronic.regions | GEX.Reads.mapped.confidently.to.transcriptome | GEX.Reads.mapped.to.genome | GEX.Reads.with.TSO | GEX.Sequenced.read.pairs | GEX.Total.genes.detected | GEX.Valid.UMIs | GEX.Valid.barcodes | Linked.genes | Linked.peaks | library_name | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1 | qmzb59ew_t11l8dzm | cellranger-arc-1.0.0 | GRCh38 | 5530 | 0.90964 | 0.05144 | 0.90915 | 0.44015 | 0.77105 | 0.75653 | 57807.81 | 8271 | 0.02425 | 132440 | 0.70263 | 0.91053 | 0.96303 | 0.96083 | NA | 319677205 | 11.26773 | 0.01051 | 0.97903 | 151554 | 0.90220 | 39397.02 | 4847.0 | 2142 | 0.68264 | 0.96989 | 0.97308 | 0.95016 | 0.22107 | 0.37121 | 0.91608 | 0.05495 | 0.48992 | 0.63371 | 0.96232 | 0.05575 | 217865541 | 28636 | 0.99949 | 0.94706 | 11028 | 45856 | BCLL_9_T_2 |
| 2 | co7dzuup_xuczw9vc | cellranger-arc-1.0.0 | GRCh38 | 5572 | 0.90972 | 0.04413 | 0.85491 | 0.42580 | 0.74570 | 0.72928 | 55415.02 | 8776 | 0.02305 | 115988 | 0.66436 | 0.89602 | 0.95478 | 0.95245 | NA | 308772502 | 10.55462 | 0.00945 | 0.97738 | 147451 | 0.91009 | 35104.20 | 4322.0 | 1997 | 0.68069 | 0.97089 | 0.97399 | 0.94723 | 0.22482 | 0.36198 | 0.92291 | 0.05291 | 0.50802 | 0.63881 | 0.96428 | 0.04315 | 195600575 | 28243 | 0.99947 | 0.94793 | 10657 | 42827 | BCLL_9_T_1 |
| 3 | wdp0p728_jf6w68km | cellranger-arc-1.0.0 | GRCh38 | 7861 | 0.90808 | 0.03255 | 0.76271 | 0.40009 | 0.68396 | 0.66529 | 31841.03 | 6958 | 0.00803 | 98142 | 0.49401 | 0.91650 | 0.96424 | 0.96264 | NA | 250302362 | 9.61085 | 0.00692 | 0.98580 | 157487 | 0.92367 | 24671.96 | 4557.0 | 2117 | 0.52518 | 0.97029 | 0.97346 | 0.95037 | 0.22588 | 0.38523 | 0.92960 | 0.05497 | 0.48940 | 0.64239 | 0.96983 | 0.05795 | 193946270 | 28847 | 0.99949 | 0.94805 | 11387 | 38686 | BCLL_8_T_2 |
| 4 | ulx1v6sz_8a2nvf1c | cellranger-arc-1.0.0 | GRCh38 | 7181 | 0.90888 | 0.03716 | 0.80033 | 0.42550 | 0.72507 | 0.70834 | 31305.98 | 7110 | 0.00855 | 108000 | 0.49722 | 0.91854 | 0.96808 | 0.96549 | NA | 224808216 | 10.26944 | 0.00778 | 0.98630 | 170819 | 0.92461 | 30537.81 | 5236.0 | 2318 | 0.55138 | 0.96930 | 0.97212 | 0.95113 | 0.22932 | 0.37099 | 0.92631 | 0.05645 | 0.49887 | 0.63400 | 0.97059 | 0.05456 | 219291981 | 29095 | 0.99950 | 0.94705 | 12085 | 42358 | BCLL_8_T_1 |
| 5 | admae8w2_89i88tvv | cellranger-arc-1.0.0 | GRCh38 | 5788 | 0.89701 | 0.02220 | 0.60932 | 0.45470 | 0.64987 | 0.62956 | 39558.28 | 7394 | 0.00290 | 74187 | 0.42259 | 0.90235 | 0.95770 | 0.95480 | NA | 228963303 | 10.45955 | 0.01080 | 0.98068 | 139918 | 0.87102 | 24146.67 | 3476.5 | 1713 | 0.56009 | 0.97060 | 0.97394 | 0.95230 | 0.17797 | 0.43160 | 0.91885 | 0.06388 | 0.42336 | 0.67030 | 0.97442 | 0.06432 | 139760911 | 27740 | 0.99953 | 0.95331 | 8430 | 30643 | BCLL_15_T_1 |
| 6 | sr20954q_yiuuoxng | cellranger-arc-1.0.0 | GRCh38 | 5845 | 0.89716 | 0.02227 | 0.61257 | 0.46364 | 0.66085 | 0.64027 | 39362.91 | 7073 | 0.00300 | 73885 | 0.42324 | 0.90302 | 0.95704 | 0.95494 | NA | 230076230 | 10.55038 | 0.01108 | 0.98122 | 140366 | 0.88381 | 30731.03 | 4005.0 | 1895 | 0.60792 | 0.96944 | 0.97246 | 0.95280 | 0.17746 | 0.43255 | 0.91632 | 0.06370 | 0.42007 | 0.66831 | 0.97285 | 0.06891 | 179622850 | 28285 | 0.99953 | 0.95371 | 8855 | 30571 | BCLL_15_T_2 |
| 7 | bs2e7lr7_mdfwypvz | cellranger-arc-1.0.0 | GRCh38 | 7460 | 0.88155 | 0.01497 | 0.44372 | 0.41225 | 0.52801 | 0.50781 | 26579.67 | 4298 | 0.01152 | 54687 | 0.26634 | 0.89865 | 0.95679 | 0.95465 | NA | 198284372 | 9.51567 | 0.00723 | 0.97758 | 19393 | 0.90605 | 21134.62 | 3354.0 | 1586 | 0.55768 | 0.96849 | 0.97215 | 0.94678 | 0.12086 | 0.35940 | 0.91054 | 0.07914 | 0.47199 | 0.70309 | 0.97493 | 0.11828 | 157664248 | 28650 | 0.99964 | 0.96253 | 5099 | 10548 | BCLL_2_T_3 |
| 8 | pd9avu0k_kf9ft6kk | cellranger-arc-1.0.0 | GRCh38 | 6919 | 0.91678 | 0.05087 | 0.91946 | 0.44305 | 0.73416 | 0.71805 | 24739.50 | 7714 | 0.00444 | 123616 | 0.44733 | 0.88822 | 0.94466 | 0.94412 | NA | 171172604 | 10.70481 | 0.01136 | 0.98087 | 144699 | 0.91809 | 28372.47 | 3828.0 | 1736 | 0.67548 | 0.96954 | 0.97278 | 0.95094 | 0.17905 | 0.41910 | 0.91917 | 0.06488 | 0.43519 | 0.66884 | 0.96896 | 0.06737 | 196309095 | 28051 | 0.99954 | 0.95416 | 9922 | 45191 | BCLL_14_T_1 |
| 9 | vuuqir4h_wfkyb5v8 | cellranger-arc-1.0.0 | GRCh38 | 6103 | 0.91469 | 0.04771 | 0.90557 | 0.43784 | 0.72039 | 0.70473 | 25950.80 | 7821 | 0.00411 | 120001 | 0.42968 | 0.89776 | 0.95826 | 0.95511 | NA | 158377744 | 10.74407 | 0.01220 | 0.97978 | 123546 | 0.90995 | 23631.59 | 3561.0 | 1659 | 0.62458 | 0.97037 | 0.97379 | 0.94388 | 0.17702 | 0.41979 | 0.92120 | 0.06485 | 0.43656 | 0.67284 | 0.96767 | 0.08440 | 144223620 | 27459 | 0.99958 | 0.95461 | 8947 | 40870 | BCLL_14_T_2 |
| 10 | kmbfo1ab_ie02b4ny | cellranger-arc-1.0.0 | GRCh38 | 10837 | 0.87842 | 0.01395 | 0.46442 | 0.36849 | 0.45577 | 0.43771 | 23010.08 | 3623 | 0.00997 | 51351 | 0.25441 | 0.90147 | 0.95567 | 0.95404 | NA | 249360235 | 8.97706 | 0.00718 | 0.98139 | 15409 | 0.92051 | 14295.57 | 2425.0 | 1291 | 0.49744 | 0.96769 | 0.97084 | 0.94821 | 0.12348 | 0.35548 | 0.90051 | 0.07662 | 0.46840 | 0.69285 | 0.97212 | 0.11530 | 154921110 | 29143 | 0.99959 | 0.96211 | 4788 | 9182 | BCLL_2_T_1 |
| 11 | ryh4el3i_biv0w7ca | cellranger-arc-1.0.0 | GRCh38 | 7910 | 0.88602 | 0.01503 | 0.45653 | 0.40099 | 0.51282 | 0.49108 | 17759.52 | 3088 | 0.01001 | 53499 | 0.21629 | 0.88597 | 0.94668 | 0.94405 | NA | 140477823 | 9.16290 | 0.00793 | 0.97915 | 13469 | 0.90946 | 20475.87 | 3298.5 | 1584 | 0.54781 | 0.97118 | 0.97343 | 0.94894 | 0.11372 | 0.35705 | 0.89950 | 0.07922 | 0.46323 | 0.69910 | 0.96995 | 0.11509 | 161964102 | 28869 | 0.99957 | 0.96247 | 4640 | 8671 | BCLL_2_T_2 |
Let us start by plotting the estimated number of cells per library. Note this this number will differ a lot from the final number of cells after applying future QC filters.
barplot_data(melt(qc_samples[,c("library_name","Estimated.number.of.cells")])) + coord_flip()
print(paste("Estimated number of cells for all the samples:", sum(qc_samples$Estimated.number.of.cells)))
## [1] "Estimated number of cells for all the samples: 77006"
Here, we can visualize the number of the features (peaks and genes) and feature linkages detected. Note that the feature linkage is calculated taking in account genes with its proximal peaks and between pairs of proximal peaks across the genome (based on correlation).
We can see how three samples BCL_2_T_3, BCL_2_T_1, BCL_2_T_2 have a low number of detected features, maybe caused by a low number of recovered nuclei These samples could be classified a priori as low quality samples.
barplot_data(melt(qc_samples[,c("library_name",
"Linked.genes",
"Linked.peaks")]))
linkage_melt <- melt(qc_samples[,c("library_name",
"Feature.linkages.detected")])
barplot_data(multiome_melted = linkage_melt)
Our first objective is to evaluated the quality and the quantity of our sequenced libraries prior to mapping.
Here, you can see the total number of sequenced read pairs assigned to the Multiome ATAC library. The optimal number of read pairs are 25,000 per cell taking in account the technical note.
qc_samples$ATAC.Mean.raw.read.pairs.per.cell <- round(qc_samples$ATAC.Mean.raw.read.pairs.per.cell,2)
ggbarplot(qc_samples,
x = "library_name",
y = "ATAC.Mean.raw.read.pairs.per.cell",
label = TRUE,
fill = "gray70",
ggtheme = theme_pubr(x.text.angle = 45,legend = "none")) +
geom_hline(yintercept = 25000, linetype='dotted', col = 'black')
Here, you can see the quality (leveraging the “Q30” variables) of sequenced read pairs assigned to the Multiome ATAC library.
multiome_melted <- melt(qc_samples[,c("library_name", "ATAC.Q30.bases.in.barcode",
"ATAC.Q30.bases.in.read.1",
"ATAC.Q30.bases.in.read.2")])
ggplot(multiome_melted, aes(variable, value)) +
geom_boxplot(fill = NA, outlier.shape = NA) +
geom_jitter(width = 0.1, height = 0) +
geom_text_repel(data = multiome_melted, aes(label = library_name, fill = "black")) +
theme_minimal() +
theme(legend.position = "none",
axis.title = element_text(size = 10),
axis.text = element_text(size = 10, color = "black"))
Fraction of read pairs with barcodes that match the whitelist after error correction. A value higher than 85% represent a high quality library.
ggbarplot(qc_samples,
x = "library_name",
y = "ATAC.Valid.barcodes",
label = TRUE,
fill = "gray70",
ggtheme = theme_pubr(x.text.angle = 45,legend = "none")) +
geom_hline(yintercept = 0.85, linetype='dotted', col = 'black')
The percentatge of the duplicated data is correlated with the library size and the sequencing depth.
correlation_plot(qc_samples,"ATAC.Percent.duplicates","ATAC.Sequenced.read.pairs")
correlation_plot(qc_samples,"ATAC.Percent.duplicates","ATAC.Number.of.peaks")
High-Quality Fragment: A read-pair with mapping quality > 30, that is not chimerically mapped, has a valid 10x barcode, and maps to any nuclear contig (not mitochondria) that contains at least one gene.
qc_atac_read_pairs <- melt(qc_samples[,c("library_name","ATAC.Mean.raw.read.pairs.per.cell", "ATAC.Median.high.quality.fragments.per.cell")])
barplot_data(qc_atac_read_pairs)
qc_atac_fraction <- melt(qc_samples[,c("library_name","ATAC.Fraction.of.high.quality.fragments.in.cells","ATAC.Fraction.of.transposition.events.in.peaks.in.cells")])
barplot_data(qc_atac_fraction)
The targetting metrics can be summarized by these 4 main score:
This number presents a high correlation with the sequencing depth specifically with the high quality fragments.
ggbarplot(qc_samples, "library_name", "ATAC.Number.of.peaks",
label = TRUE,
fill = "gray70",
ggtheme = theme_pubr(x.text.angle = 45,legend = "none"),
position = position_dodge(0.9))
correlation_plot(qc_samples,"ATAC.Number.of.peaks",
"ATAC.Median.high.quality.fragments.per.cell")
The fraction of high quality fragments in cell are expected to be higher 40%. The fraction of transposition events that fall within peaks > 25%
multiome_melted <- melt(qc_samples[,c("library_name", "ATAC.Fraction.of.high.quality.fragments.in.cells",
"ATAC.Fraction.of.high.quality.fragments.overlapping.peaks", "ATAC.Fraction.of.high.quality.fragments.overlapping.TSS")])
multiome_melted$value = round(multiome_melted$value,2)
ggbarplot(multiome_melted, "library_name", "value",
fill = "variable", color = "variable", palette = "Paired",
label = TRUE,
ggtheme = theme_pubr(x.text.angle = 45),
position = position_dodge(0.9))
This fraction is quite low in all samples. We do not have a reference value to be able to compare them. However, the examples we reviewed had a value of around 2%.
ggplot(melt(qc_samples[,c("library_name",
"ATAC.Fraction.of.genome.in.peaks")]), aes(variable, value)) +
geom_boxplot(fill = NA, outlier.shape = NA) +
geom_jitter(width = 0.1, height = 0) +
theme_minimal() +
theme(legend.position = "none",
axis.title = element_text(size = 10),
axis.text = element_text(size = 10, color = "black"))
It is expected to obtain a large enrichment around TSS (> 5% EV, red dashed line), as these regions are known to show a high degree of chromatin accessibility compared to the flanking regions.
data_tss <- qc_samples[, c("library_name", "ATAC.TSS.enrichment.score")]
data_tss.melt <- melt(data_tss)
barplot_data(data_tss.melt) + geom_hline(yintercept = 5, linetype = 2, color = "red")
median_genes_gg <- melt(qc_samples[,c("library_name","GEX.Median.genes.per.cell")])
barplot_data(median_genes_gg)
Here, you can see the total number of sequenced read pairs assigned to the Multiome Gene Expression library. The optimal number of read pairs are ~20,000 per cell taking in account the technical note.
qc_samples$GEX.Mean.raw.reads.per.cell <- round(qc_samples$GEX.Mean.raw.reads.per.cell,2)
ggbarplot(qc_samples,
x = "library_name",
y = "GEX.Mean.raw.reads.per.cell",
label = TRUE,
fill = "gray70",
ggtheme = theme_pubr(x.text.angle = 45)) +
geom_hline(yintercept = 20000, linetype='dotted', col = 'black')
Here, you can see the quality (leveraging the “Q30” variables) of sequenced read pairs assigned to the Multiome Gene Expresion library.
multiome_melted <- melt(qc_samples[,c("library_name", "GEX.Q30.bases.in.barcode",
"GEX.Q30.bases.in.UMI",
"GEX.Q30.bases.in.read.2")])
ggplot(multiome_melted, aes(variable, value)) +
geom_boxplot(fill = NA, outlier.shape = NA) +
geom_jitter(width = 0.1, height = 0) +
geom_text_repel(data = multiome_melted, aes(label = library_name, fill = "black")) +
theme_minimal() +
theme(legend.position = "none",
axis.title = element_text(size = 10),
axis.text = element_text(size = 10, color = "black"))
Fraction of read pairs with barcodes that match the whitelist after error correction. A value higher than 80% represent a high quality library,
ggbarplot(qc_samples,
x = "library_name",
y = "GEX.Valid.barcodes",
label = TRUE,
fill = "gray70",
ggtheme = theme_pubr(x.text.angle = 45,legend = "none")) +
geom_hline(yintercept = 0.85, linetype='dotted', col = 'black')
The percentatge of the duplicated data is correlated with the library size and the sequencing depth.
correlation_plot(qc_samples,"GEX.Percent.duplicates","GEX.Sequenced.read.pairs")
correlation_plot(qc_samples,"GEX.Percent.duplicates","GEX.Total.genes.detected")
Secondly, we will assess the quality of cellranger’s mapping by comparing the percentage of reads mapping to the genome, intergenic regions, intronic and exonic regions across libraries. Reads mapping to intergenic regions suggest contamination of ambient DNA, while reads mapping to intronic regions may come from pre-mRNAs or mature splice isoforms that retain the introns.
The fraction of sequenced reads that map to a unique gene in the transcriptome is expected to be higher than the 50%.
qc_samples$GEX.Fraction.of.transcriptomic.reads.in.cells <- round(qc_samples$GEX.Fraction.of.transcriptomic.reads.in.cells,2)
ggbarplot(qc_samples,
x = "library_name",
y = "GEX.Fraction.of.transcriptomic.reads.in.cells",
fill = "gray70",
label = TRUE,
ggtheme = theme_pubr(x.text.angle = 45)) +
geom_hline(yintercept = 0.50, linetype='dotted', col = 'black')
multiome_melted <- melt(qc_samples[,c("library_name", "GEX.Reads.mapped.to.genome", "GEX.Reads.mapped.confidently.to.genome")])
multiome_melted$value <- round(multiome_melted$value,2)
ggbarplot(multiome_melted,
x = "library_name",
y = "value",
fill = "variable", palette = "Paired", label = TRUE,
ggtheme = theme_pubr(x.text.angle = 45),
position = position_dodge(0.9))
The fraction of sequenced reads that map to a unique gene in the transcriptome is expected to be higher than 50% (blue dashed line). However, the percentatge of reads mapped confidently to intergenic regions and reads mapped antisense to gene are expected to be lower than 30% (red dashed line).
mapping_qc_vars <- c(
"library_name",
"GEX.Reads.mapped.confidently.to.intergenic.regions",
"GEX.Reads.mapped.confidently.to.intronic.regions",
"GEX.Reads.mapped.confidently.to.exonic.regions",
"GEX.Reads.mapped.confidently.to.transcriptome",
"GEX.Reads.mapped.antisense.to.gene")
mapping_qc_gg.melt <- melt(qc_samples[,mapping_qc_vars])
ggplot(mapping_qc_gg.melt, aes(variable, value)) +
geom_boxplot(fill = NA, outlier.shape = NA) +
geom_jitter(width = 0.1, height = 0) +
theme_minimal() +
theme(legend.position = "none",
axis.title = element_text(size = 12),
axis.text = element_text(size = 12, color = "black")) + coord_flip() +
geom_hline(yintercept = c(0.3,0.5), col = c("red","blue"), linetype='dotted')
Thirdly, we will plot the number of detected genes per library as a function of the total reads sequenced. We know that this function reaches a plateau, whereby more sequenced reads does not result in more detected genes. In those scenarios, we say that we have sequenced until saturation:
ggplot(qc_samples, aes(GEX.Sequenced.read.pairs, GEX.Total.genes.detected)) +
geom_point() +
geom_text_repel(data = qc_samples, aes(label = library_name), color = "black") +
labs(x = "Number of Reads", y = "Total Genes Detected", color = "") +theme_classic() +
theme(axis.title = element_text(size = 13, color = "black"),
axis.text = element_text(size = 12, color = "black"),
legend.text = element_text(size = 12, color = "black"))
sessionInfo()
## R version 4.0.3 (2020-10-10)
## Platform: x86_64-apple-darwin13.4.0 (64-bit)
## Running under: macOS Big Sur 10.16
##
## Matrix products: default
## BLAS/LAPACK: /Users/pauli/opt/anaconda3/envs/Tonsil_atlas/lib/libopenblasp-r0.3.10.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] plotly_4.9.2.1 ggrepel_0.8.2 kableExtra_1.3.1 knitr_1.30 data.table_1.13.2 reshape2_1.4.4 plyr_1.8.6 ggpubr_0.4.0 ggplot2_3.3.2 BiocStyle_2.16.1
##
## loaded via a namespace (and not attached):
## [1] Rcpp_1.0.5 lattice_0.20-41 here_1.0.1 tidyr_1.1.2 rprojroot_2.0.2 digest_0.6.27 R6_2.5.0 cellranger_1.1.0 backports_1.2.0 evaluate_0.14 httr_1.4.2 highr_0.8 pillar_1.4.6 rlang_0.4.8 lazyeval_0.2.2 curl_4.3 readxl_1.3.1 rstudioapi_0.12 car_3.0-10 Matrix_1.2-18 rmarkdown_2.5 splines_4.0.3 labeling_0.4.2 webshot_0.5.2 stringr_1.4.0 foreign_0.8-80 htmlwidgets_1.5.2 munsell_0.5.0 broom_0.7.2 compiler_4.0.3 xfun_0.18 pkgconfig_2.0.3 mgcv_1.8-33 htmltools_0.5.0 tidyselect_1.1.0 tibble_3.0.4 bookdown_0.21 rio_0.5.16 viridisLite_0.3.0 crayon_1.3.4 dplyr_1.0.2 withr_2.3.0 grid_4.0.3 nlme_3.1-150 jsonlite_1.7.1 gtable_0.3.0 lifecycle_0.2.0 magrittr_1.5 scales_1.1.1 zip_2.1.1 stringi_1.5.3 carData_3.0-4 farver_2.0.3 ggsignif_0.6.0 xml2_1.3.2 ellipsis_0.3.1 generics_0.1.0 vctrs_0.3.4 openxlsx_4.2.3
## [60] RColorBrewer_1.1-2 tools_4.0.3 forcats_0.5.0 glue_1.4.2 purrr_0.3.4 hms_0.5.3 abind_1.4-5 yaml_2.2.1 colorspace_2.0-0 BiocManager_1.30.10 rstatix_0.6.0 rvest_0.3.6 haven_2.3.1